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1.  Introduction 

In  order  to  increase  communication  between  investigators  who  are  modeling  upper 
ocean  biological/pbysical  processes,  a  workshop  was  convened  in  Woods  Hole  on  June  7- 
12,  1993  (see  Appendix  I  for  names  and  addresses  of  participants).  This  workshop  was 
part  of  our  on-going  URIP  project  entitled  “Modehng  Biological-Physical  Interactions;  A 
Population  Biological  Approach”  sponsored  by  ONR  (Grant  N00014-92-J-1527).  The  two 
principal  goals  of  the  workshop  were  to:  1)  identify  critical  problems  related  to  mixed-layer 
biological-physical  models,  and  2)  develop  approaches  for  solving  these  problems. 

The  workshop  was  organized  into  two  parts  to  address  these  goals.  The  first  part, 
held  over  the  first  day  and  a  half,  included  three  overview  presentations  given  in  plenary 
followed  by  working  groups,  organized  along  disciphnary  lines,  to  identify  critical  issues. 
The  second  part  of  the  workshop  consisted  of  working  groups,  organized  across  discipHnes, 
using  “hands-on”  modehng  to  address  critical  aspects  of  coupled  biological-physical  models. 

Overview  presentations  were  given  during  the  first  morning  to  provide  background 
information  on  physical  mixed-layer  models  (Kenneth  Denman),  biological  models  with 
age-structured  herbivore  populations  (John  Steele),  and  microbial-loop  dynamics  (David 
Caron).  Following  these  presentations,  the  working  groups  were  organized  to  address 
problems  associated  with  1)  physical  modehng,  2)  structured  herbivore  vs  bulk  nutrient- 
phytoplankton-zooplankton  models,  3)  food-web  models,  and  4)  couphng  complex  biology 
to  2-  and  3-D  physical  models. 

After  plenary  reports  of  these  first-format  working  groups,  the  second-format  working 
groups  were  formed  to  begin  modehng  in  the  areas  of  1)  couphng  simple  food- web  models  to 
several  different  physical  models,  2)  development  of  NPZ  models  which  include  microbial- 
loop  dynamics,  and  3)  development  of  structured-zooplankton  models  and  incorporation 
into  NPZ  models. 

An  annotated  bibhography  of  existing  mixed-layer  models  was  provided  by  David 
Archer  and  is  given  in  Appendix  11.  Brief  descriptions  of  three  coupled  biological/physical 
mixed  layer  models  is  given  in  Appendix  III. 

2.  General  Reviews-Critical  Issues  (from  Ist-format  Working  Groups) 

2.1  Physics  of  Mixed  Layers 

Archer,  Chen,  Denman,  Doney,  Gawarkiewcz  (Rapporteur),  Glover  (Chair),  Hood 

•  The  biological  questions  which  are  being  addressed  should  affect  the  choice  of  model 
being  used. 

•  Models  have  clearly  defined  weaknesses  which  must  be  taken  into  consideration  before 
use.  For  example,  P-W-P  assumes  no  vertical  velocity  shear  is  present  in  the  mixed 
layer,  while  observations  show  such  shear  is  frequently  present. 

•  We  should  be  thinking  about  upper  ocean  models  in  general  (i.e.  including  the  ther- 
mochne)  rather  than  strictly  hmiting  ourselves  to  calculating  the  ’’mixed  layer  depth”. 
This  is  particularly  important  when  considering  heat  flux  divergences  in  1-D  models 
(what  depth  do  you  ’’advect”  the  heat  away  when  modelling  an  annual  cycle?).  Lat¬ 
eral  mixing  processes  may  be  important  in  the  thermocline  in  distributing  nutrients 
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along  isopycnals.  The  vertical  transport  of  nutrients  into  the  base  of  the  mixed  layer 
needs  to  be  considered. 

•  Observations  of  Jenkins  clearly  show  an  enhanced  vertical  flux  of  nutrients  into  the 
euphotic  zones  which  current  models  using  standard  parameterizations  cannot  explain. 
Is  this  due  to  mesoscale  variabihty  and  eddy  pumping?  What  should  be  done  with 
1-D  parameterizations  if  mesoscale  effects  are  important? 

•  What  role  do  large  and  small-scale  fronts  play  in  understanding  mixed  layer  dynamics? 
How  important  are  horizontal  inhomogeneities  in  understanding  the  performance  of  1- 
D  models  when  compared  with  time  series  of  observations?  How  should  mixing  within 
frontal  zones  be  parameterized?  More  observations  with  good  horizontal  resolution 
wiU  be  necessary  to  make  progress  in  this  area. 

•  What  level  of  complexity  is  needed  in  prescribing  optical  fields?  Under  what  conditions 
do  the  absorption  characteristics  affect  the  temperature  distributions  in  the  vertical? 

•  There  are  major  questions  regarding  the  Lagrangian  descriptions  of  flow  fields  in  the 
upper  ocean.  In  particular,  the  Lagrangian  behavior  of  particles  in  fields  where  the 
mixing  coefficients  are  changing  with  depth  is  not  well  understood.  Good  Lagrangian 
information  is  necessary  to  describe  light  histories  of  particles  as  well  as  predator-prey 
contact  rates. 

•  How  will  the  use  of  Large  Eddy  Simulations  affect  our  understanding  of  the  mixed 
layer?  These  models  seem  to  give  results  that  are  very  different  from  both  P-W-P  as 
well  as  MeUor-Yamada  Level  2.5. 

•  The  integration  of  data  from  various  acoustic  observational  techniques  which  resolve 
the  turbulent  eddies  in  the  mixed  layer  should  lead  to  much  better  parameterizations 
of  the  characteristics  of  the  turbulence.  Hopefully,  this  will  occur  over  the  next  5  years 
or  so. 

•  In  the  coastal  regions,  the  interaction  of  the  bottom  and  surface  boundary  layers  is 
not  weU  understood.  Much  more  work  in  the  future  wiU  be  necessary  to  resolve  this 
issue. 

•  Surface  forcing  in  regions  in  which  no  buoy  data  is  available  is  a  major  problem. 
Wind  products  tend  to  severely  underestimate  maximum  wind  velocities,  which  leads 
to  major  differences  in  surface  fluxes. 

•  Processes  associated  with  the  surface  wave  field  are  not  included  in  most  present 
day  mixed  layer  models.  The  effects  of  Langmuir  cells  and  wave  dissipation  on  the 
turbulence  fields  will  have  to  be  addressed  in  the  future. 

Conclusions 

•  Horizontal  variability  is  not  weU  known. 

•  Lagrangian  nature  of  flows  are  not  well  known. 

•  Acoustic  observational  techniques  will  help  resolve  turbulence  fields. 

2.2  NPZ  vs  NP-( Structured- Z) 

Armstrong,  BoUens,  CasweU  (Chair),  Frost,  Lewis  (Rapporteur),  Steele 

•  Why  structure  only  herbivore  populations? 
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-  No  reason,  in  fact  the  group  decided  in  many  cases  that  structure  could  be  re¬ 
quired  in  a  number  of  different  trophic  levels. 

What  is  the  purpose  of  structure? 

-  This  was  one  of  the  primary  questions  that  came  around,  as  it  is  possible  to  add 
structure  to  a  wide  variety  of  popiilations  and  trophic  levels  while  accomphsh- 
ing  nothing  more  than  tying  up  your  computer.  Apphcations  in  which  stage- 
structured  models  may  be  needed  include: 

-  Fisheries  models  in  which  it  is  necessary  to  know  the  abundance  of  a  wide  variety 
of  zooplankton  species;  prey  for  the  larval  fish. 

-  Data  comparison  and  validation.  However,  the  converse  of  this  is  that  there  is 
little  point  in  running  a  stage  structured  model  where  you  have  no  information 
concerning  the  population  structure. 

-  For  comparison  to  NPZ  models  to  determine  whether  the  introduction  of  stage 
structure  has  a  significant  effect  on  the  dynamics. 

-  To  introduce  stage  specific  behavior. 

-  To  represent  the  microbial  loop  and  other  multispecies  models. 

-  Time  scales  of  response  by  a  species  with  a  long  growth  period  may  require  stage 
structure  to  be  adequately  represented. 

-  Abundances  may  respond  with  a  periodicity  determined  by  the  generation  times. 

Problems  with  such  stage  models  include: 

-  You  may  need  multiple  phytoplankton  populations  to  support  a  multi-  species 

grazing  population. 

-  Truncation  of  the  trophic  models  to  the  highest  trophic  level  of  relevance. 

Vertical  structure: 

-  Need  the  physical  dynamics  to  transport  plankton,  but  you  may  not  need  a  full 
one-D  biological  model  to  couple  to  the  physics. 

-  Behavioral  effects. 

Incorporation  of  physical  models. 

-  May  not  need  full  resolution  of  PWP  model. 

-  Should  explore  sensitivity  of  models  to  increased  physical  complexity. 

Suggestions:  _ 

-  Use  output  of  physical  model  to  drive  biology:  use  temperature  and  diffusivity 

values,  time-smoothed 

-  Time  smooth  and  space  smooth  spatial  model  at  varying  scales  to  determine  what 
effect  the  reduced  complexity  will  have  on  the  biology. 

-  Use  100  point  PWP  results  on  100,50,20,  etc.  level  biological  model. 

-  Use  time  averaged  PWP  fields,  determine  effect  of  high  frequency  activity. 

-  Combine  the  organisms  of  Caron’s  work  with  the  model  written  by  Armstrong. 

Food  web  models  and  validation  with  data 

Caron,  Landry  (Chair),  Landsteiner  (Rapporteur),  Moisan,  Sarmiento 

Critical  Task: 

-  Implementation  of  food- web  models  with  appropriate  microbial  community  struc¬ 
ture  and  interactions. 


•  Minimal  structure  of  a  microbial  loop  model: 

-  Three  phytoplankton  (pico,  nano,  and  micro) 

-  Heterotrophic  bacteria 

-  Two  levels  of  protozoan  consumers  (nano  and  micro) 

-  DOM  flows  from  phytoplankton  excretion,  sloppy  feeding,  decay  of  protozoan 
feces,  and  mortality. 

•  Relevance  of  microbial  dynamics  in  food  webs:  -  Basic  structure  of  lower  trophic  levels 
in  all  upper  ocean  environments. 

—  Allows  partitioning  of  primary  production  among  size  classes,  a  major  distin¬ 
guishing  characteristic  of  different  environments  or  seasons. 

-  Essential  to  understanding  the  basic  mechanisms  of  phytoplankton  control  (i.e., 
grazing  vs  nutrient  hmitation). 

-  Provides  a  realistic  resource  environment  for  modehng  of  macro-consumers. 

-  Provides  a  versatile  format  for  other  desirable  model  outputs: 

-  new  production  ratio  (recycling  pathway) 

-  multi-elemental  cycling 

-  vertical  fluxes  (sinking  vs.  non-sinking  feces) 

•  Necessary  and  existing  data: 

-  standing  stocks  and  elemental  composition  of  all  populations 

-  size-fractionated  primary  production 

-  protozoan  grazing  rates  and  behaviors 

-  protozoan  growth  efficiencies 

•  Problems: 

-  Thresholds 

-  Closure  through  metazoans  (multiple  species  and  stage-structure) 

-  AUometric  relationships 

-  Depth  structure  (controlled  through  fight  and  nutrient  kinetics) 

-  Selective  feeding  on  multiple-size  classes 

2.4  Problems  in  use  of  1-D  coupled  biological/physical  models 

Davis,  Ffierl  (Chair),  Franks  (Rapporteur),  Levin,  McGilficuddy,  Olson 

•  Critical  Issues 

—  Horizontal  advection  necessary  over  long  time-scales  to  balance  deep  heat  flux, 
nutrient  flux 

_  Timesc2iles  of  mixing  in  the  model  vs.  biological  rates  -  how  to  parameterize 
faster  rates  (e.g.  photoadaptation)  in  mixed  layers 

-  How  to  parameterize  patchiness  and  variability  of  biological  fields 

-  How  does  variance  propagate  trophicaUy,  temporally  -  can  we  use  statistics  of 
biological  distributions  as  variables  in  the  model? 

-  When  is  use  of  PWP  model  appropriate? 

-  short-term  variation  in  seasonal  thermocfine;  days- weeks 

-  if  coupled  to  ’’good”  GCM,  may  be  able  to  extend  its  usefulness 
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—  need  to  consider  importance  of  remnant  layers  to  vertical  fluxes  over  long 
times  scales.  May  necessitate  forced  boundary  conditions,  or  coupling  to 
2-D  or  3-D  PE  model. 

•  Suggestions: 

-  use  a  number  of  different  models  of  the  same  physical/  biological  problem  (sensu 
Martin)  to  assess  utility,  and  observe  differences. 

-  use  large-eddy  simulation  models  for  more  realistic  2-D  or  3-D  representation  of 

smaller  scales 

-  use  Deardorff-type  model  (4th  order  closure)  rather  than  PWP 

3.  Reports  of  Work  (2nd  format  Working  Groups) 

3.1  Comparing  mixed  layer  models  containing  simple  plankton  dynamics 

Chen,  Denman  (Chair),  Doney,  Fherl,  Franks,  Gawarkiewcz,  Glover,  Hood 

3.1.1  Goal 

To  compare  the  behavior  of  several  models  of  mixed  layer  physics  coupled  to  simple 

Nutrient-Phytoplankton-Zooplankton-Detritus  (NPZD)  models. 

3.1.2  Physical  Mixed  layer  models 

The  group  concentrated  primarily  on  completing  the  comparison  and  couphng  of 
mixed  layer  models  to  a  simple  NPZ  model,  all  set  up  beforehand  by  Glenn  Fherl.  In 
addition  S.  Doney  and  D.  Glover  set  up  a  mixed  layer  model  being  developed  at  NCAR 
with  an  NPZD  model,  and  R.  Hood  explored  the  behavior  of  an  NPZD  model  coupled  to 
CONV  and  PWP  (described  below).  The  following  mixed  layer  models  were  included  in 
the  comparison: 

•  CONV  -  This  model  was  the  simplest,  consisting  of  a  mixed  layer  depth  of  10  m  or 
else  one  determined  by  simple  convective  adjustment,  whichever  was  deeper.  Each 
day  at  noon,  the  mid  retreated  to  10  m  and  each  night  it  reached  a  depth  ranging 
from  about  20  m  in  summer  to  80  m  for  1-2  months  in  the  winter  (which  was  the 
bottom  of  the  model).  The  annual  cycle  was  more  reahstic  than  in  the  plots  to  follow, 
which  track  the  mid  at  midnight,  several  hours  before  the  deepest  mixing  each  night. 

•  PWP  (Price- WeUer-Pinkel)  -  This  model  consists  of  a  mixed  layer  where  the  buoyancy 
is  mixed  completely  each  time  step  (Price,  Weller  and  Pinkel,  1986).  The  thickness 
of  the  mixed  layer  is  determined  by  a  bulk  Richardson  number  closure  (Rb  ^  0.65), 
followed  by  a  smoothing  below  the  mixed  layer  to  reheve  any  shear  instabihties  such 
that  the  local  gradient  Richardson  number  Rg  >  0.25.  This  version,  coded  by  Jim 
Price  and  modified  by  Glenn  Flierl,  also  has  an  option  for  a  background  diffusion 
which  was  not  used  except  by  R.  Hood  who  was  using  a  version  modified  by  Hood 
and  Olson  to  correct  problems  with  the  coding  of  the  background  diffusion. 

•  MY  (MeUor-Yamada)  -  THs  model  is  based  on  a  level-2  turbulent  closure  scheme 
outhned  in  MeUor  and  Yamada  (1974)  and  first  implemented  by  MeUor  and  Durbin 
(1975).  The  vertical  mixing  at  any  depth  is  determined  by  a  turbulent  eddy  coefficient 
that  is  a  function  of  the  mean  turbulent  kinetic  energy  (TKE),  a  length  scale  taken  to 
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be  the  first  moment  of  the  TKE,  and  a  stabihty  coefficient  that  is  calculated  according 
to  some  laboratory-determined  empirical  equations  that  are  functions  of  the  local  flux 
Richardson  number  Rf.  Thus,  an  output  of  the  model  is  a  vertical  profile  of  the 
coefficient  of  vertical  turbulent  diffusion  of  heat,  Kriz),  which  is  the  key  output  for 
coupHng  with  a  biological  model.  This  version,  coded  by  Patrice  Klein  and  modified 
by  Jim  Price,  also  has  a  constant  background  diffusion. 

These  three  models  were  set  up  for  the  workshop  by  Glenn  Flierl,  with  a  common 
module  for  wind  and  heat  forcing  and  a  common  biological  NPZ  module  for  couphng.  A 
simple-to-use  set  of  plotting  commands  was  also  set  up  so  that  comparison  on  common 
axes  at  the  same  size  could  easily  by  made. 

•  NCAR  -  This  model  was  implemented  and  coupled  to  an  NPZD  at  the  meeting  by 
S.  Doney  and  D.  Glover.  This  is  a  one-dimensional  coupled  biological-physical  model 
built  at  NCAR  as  the  first  part  of  a  longer  term  project  for  developing  a  3-D,  global 
coupled  model.  The  physical  model,  which  is  based  on  models  of  the  atmospheric 
planetary  boundary  layer,  differs  from  existing  ocean  models,  and  a  complete  descrip¬ 
tion  can  be  found  in  Large  et  al.  (1993).  The  model  is  a  non-local  parameterization  of 
the  oceanic  boundary  layer  based  on  the  concept  that  the  profiles  of  eddy  diffusivity 
in  the  planetary  boundary  layer  (either  atmosphere  or  ocean)  are  similar  in  shape 
if  the  appropriate  scahng  is  appHed  (boundary  layer  depth,  surface  wind  stress  and 
surface  buoyancy  flux).  The  shape  parameterization  used  in  the  model  is  based  on 
Large  Eddy  Simulation  (LES)  results  (Deardorff,  1972).  Constraints  are  also  apphed 
to  the  parameterization  such  that  the  eddy  diffusivity  in  the  surface  layer  matches  the 
results  of  similarity  theory.  The  depth  of  the  ocean  boundary  layer  depth  is  computed 

based  on  a  bulk  Richardson  number  criteria. 

The  physical  model  can  be  run  using  either  prescribed  surface  fluxes  or  atmospheric 
state  variables  (e.g.  air  temperature  and  humidity,  wind  speed);  atmospheric  forcing 
can  be  either  specified  using  simple  harmonics  or  read  in  from  a  data  file.  The  daily 
cycle  of  solar  radiation  is  computed  using  a  geometric  solar  model  and  the  cloud  model 

of  Smith  and  Dobson  (1984). 

3.1.3  Biological  models 

The  NPZ  model  set  up  on  the  computer  by  Fherl  for  couphng  to  the  mixed  layer 
models  was  a  relatively  simple  model  with  some  non-standard  terms  such  as  the  extra 
P  factor  in  the  grazing  term  (to  replace  a  threshold  value  in  the  Ivlev  formulation)  the 
Michaehs-Menten  phytoplankton  growth  dependence  on  fight,  and  the  Z  dependence 
death  term  for  the  zooplankton.  The  model  consisted  of  the  foUowing  equations: 

RL  =  V  - —P  -  Rm>^PZ[l  -  exp(-AP)]  -  bP 

Dt  Io{t)I{z)  +  13  kg  +  N 


=  'fR^XPZ[l  -  exp(-AP)]  -  dZ"^ 


DN 

Dt 


DP  DZ 


(1) 


Dt  Dt 
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where  I(z)  is  the  hght  intensity  function 


I{z)  =  exp(z/;52) 

The  surface  Hght  intensity  Io{t)  varied  with  time  of  day  and  season.  Parameters  values 
and  initial  conditions  are  given  in  Table  I. 


Table  I.  Parameter  values  used  in  the  NPZ  model  set  up  by  Glenn  FHerl 


Parameter 

Value 

Description 

I 

_ 

irradiance 

Vm 

2.0 

maximum  phytoplankton  growth  rate 

kg 

0.1 

nutrient  half  saturation  coefficent 

Rm 

0.35 

maximum  zooplankton  growth  rate 

A 

1.0 

Ivlev  coefficent 

b 

0.05 

phytoplankton  mortaHty 

7 

0.7 

egestion  fraction 

d 

0.2 

zooplankton  mortaHty 

8 

0.5-1.5 

Hght  half-saturation  coefficient 

^2 

20. 

Hght  e-folding  depth 

State  Variables: 

II 

Zo  =  0.35  Nt>  =  1.95 

Two  other  biological  models  were  used  as  weU.  The  NPZD  model  used  by  Hood 
consisted  of  a  detrital  compartment  (D)  added  to  the  FHerl  and  Davis  (1993)  NPZ  model 
(see  section  3.1.7). 

The  biological  model  used  by  Doney/Glover  was  similar  and  was  a  nitrogen  based  flow 
model  adapted  from  Fasham’s  model  of  the  noixed  layer  (Fasham  et  al.,  1993).  The  model 
is  quite  flexible  in  that  the  flow  pathways  and  model  parameters  are  specified  as  part  of 
an  input  file  at  run  time.  Currently,  the  model  is  formulated  as  NPZD  model  and  includes 
both  turbulent  transport  and  sinking  of  the  biological  components.  CoupHng  between  the 
biological  model  and  the  physical  model  via  the  absorption  coefficient  for  solar  radiation 
can  also  be  included. 

The  general  equation  for  any  biological  scalar  at  a  particular  grid  level  is. 

^  ^  _  dWsinking-X; 

dt  dz  dz 

where  is  the  turbulent  flux  of  X  and  ^sinking  is  the  sinking  velocity  of  species  X. 
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The  biological  production  and  consumption  terms  for  the  NPZD  model  are: 


dt 


Vmil) 


PN 
k,  +  N 


-Rm{l 


e-^)Z  -  eP 


^  =  jRm{l  -  e-^)Z  -  gZ 

dt 


dN 

dt 


-Vm{I) 


PN 

k,+N 


+  mD 


=  (1  _  7)i2^(l  -  e-^)Z  +  eP-mD  +  gZ 

dt 


where  parameter  values  are  given  in  Table  11. 


(2) 


Table  II.  Parameter  values  used  in  the  Doney  and  Glover  NPZD  model. 


Parameter 

Value 

Description 

I 

— 

irradiance 

Vm 

2.0 

maximum  phytoplankton  growth  rate 

kg 

0.1 

nutrient  half  saturation  coefficent 

Rm 

1.0 

maximum  zooplankton  growth  rate 

A 

0.3 

Ivlev  coefficent 

e 

0.1 

phytoplankton  mortality 

7 

0.7 

egestion  fraction 

g 

0.2 

zooplankton  mortahty 

m 

1.0 

detrital  reminerahzation  rate 

3.1.4  Tasks  attempted  and  completed 

•  The  group  explored  and  compared  the  behavior  of  the  CONV,  PWP  and  MY  mixed 
layer  models.  They  compared  the  shape  of  the  vertical  temperature  profile,  the  annual 
cycles  of  surface  temperature  and  mixed  layer  depth,  conservation  of  heat  and  the 
vertical  diffusion  characteristics. 

•  The  NPZD  models  were  coupled  first  with  CONV  and  parameters  were  adjusted  to 
find  ’stable’  behavior,  i.e.  behavior  that  did  not  include  inherently  biological  hmit 
cycles,  but  rather  annual  cycles  controlled  by  a  combination  of  the  annual  cycles  in 
solar  radiation  and  mixed  layer  behavior. 

•  For  the  ’stable’  parameter  sets,  the  NPZ  model  was  coupled  to  MY  and  the  behavior 
compared  with  that  of  the  coupled  CONV  -  NPZ  model. 
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•  Similar  runs  were  conducted  with  the  NPZ  model  coupled  to  the  PWP  model 

•  The  NCAR  model  was  implemented  with  an  NPZD  model  and  tested  with  several 

multiyear  simulations  to  demonstrate  the  influence  of  the  detritus  loop. 

•  An  NPZD  model  was  coupled  to  the  CONV  model  and  parameter  space  explored.  A 

modifled  NPZD  model  was  coupled  with  the  Hood/ Olson  version  of  PWP. 

3.1.5  Results  of  the  comparison  of  mixed  layer  models 

The  PWP  and  CONV  models  behave  very  similarly,  in  terms  of  the  annual  cycle  of 
sea  surface  temperature,  mixed  layer  depth,  and  temperature  profile. 

Initially,  the  MY  model  did  not  conserve  heat,  but  rather  it  added  significant  heat 
beyond  the  net  heat  input  across  the  air-sea  boundary,  which  was  evident  in  both  the 
annual  SST  cycle  and  the  vertical  profiles  departing  from  those  of  CONV  and  PWP.  (Figs. 

1  and  2). 

MY  includes  a  bottom  boundary  layer  and  a  background  constant  diffusion,  initially 
set  at  1.34  x  10~®  m^  s“^.  Over  an  annual  cycle  a  large  amount  of  heat  was  diffused  from 
the  mixed  layer  into  the  layer  below,  especially  the  bottom  well-mixed  layer.  The  extent 
of  this  diffusion  is  shown  in  a  2-year  run  where  the  surface  heat  exchanges  were  turned  off 
(Fig.  3).  The  background  diffusion  dominated  the  modelled  diffusion  resulting  from  the 
level-2  closure  as  shown  in  Fig.  4,  which  is  a  plot  of  the  Kt(z)  profile  on  a  logio  plot  after 

200,  360  and  850  modelled  days  (day  1  =  120). 

The  non-conservation  of  heat  was  traced  to  a  function  that  was  applied  each  timestep 
to  keep  the  odd  and  even  time  leapfrog  solutions  from  diverging  from  each  other.  G.  FHerl 
replaced  the  scheme  with  a  1/2  timestep  Euler  solution  every  100  timesteps,  i.e.  twice  a 
day.  The  annual  cycles  of  SST  in  CONV  and  MY  then  were  almost  identical  (Fig.  5), 
where  the  background  diffusion  and  the  bottom  boundary  layer  in  MY  were  both  turned 

off. 

Another  feature  of  MY  is  that  convective  mixing  due  to  heat  loss  at  the  sea  surface 
takes  several  timesteps  to  reach  full  adjustment  for  a  deep  mixed  layer,  unhke  PWP  (and 
CONV)  where  complete  mixing  occurs  each  time  step.  This  behavior  will  generate  propa¬ 
gating  waves  when  PWP  is  embedded  in  a  circulation  model  and  adjacent  gridpoints  mix 
suddenly  to  different  depths.  However,  the  MY  model  appeared  to  take  close  to  10  times 
the  computer  time  of  PWP,  so  it  may  not  be  practical  for  inclusion  in  large  high  resolution 
3D  circulation  models. 

The  NCAR  mixed  layer  model,  based  on  a  fitted  similarity  solution  for  the  vertical 
profiles  of  Kt,  is  an  attractive  option.  Its  ability  to  simulate  a  stable  annual  cycle  in  the 
mixed  layer  is  shown  in  Fig.  6. 

As  part  of  the  workshop,  the  NCAR  model  was  transferred  to  and  adapted  to  run  on 
the  NASA  EOS  interdisciplinary  team  SCF  SGI  at  WHOI.  Using  a  set  of  simple,  harmonic 
forcing  functions,  we  generated  a  physical  model  solution  for  a  Sargasso  Sea  like  case  (Fig. 
6a  and  b).  The  model  produced  deep,  convective  mixing  during  winter  and  a  shallow 
mixed  layer  and  seasonal  thermochne  of  approximately  the  correct  amplitude  during  the 
summer.  A  distinct  and  stable  seasonal  cycle  developed  after  one  year.  Notably,  the  model 
retains  a  thick  layer  of  mode  water  at  about  eighteen  degrees. 

3.1.6  Results  with  coupled  NPZD  and  mixed  layer  models 
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6.  NCAR  K-profile  model  run  for  three  years  at  a  Sargasso  Sea  like  location 
(64° W,  32°N)  showing  model  stabihty  for  both  (a)  temperature  (°C)  and  (b) 
phytoplankton  (mM  N  m“*). 
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The  NPZ  model  coded  by  G.  Flierl  was  coupled  first  to  the  CONV  model  and  parame¬ 
ter  space  was  explored.  Initially,  interesting  behavior  resulted  where  several  Lotka-Volterra 
hke  hmit  cycles  occurred  each  spring-summer-autumn  period,  exhausting  nutrients  briefiy 
each  cycle.  These  were  slowed  down,  by  changing  the  parameters,  to  about  one  per  year. 
The  spring  phytoplankton  bloom  progressed  downward  with  time  until  it  got  below  the 
mixed  layer  where  for  some  parameter  sets,  hmit  cycles  occurred  as  the  peak  in  P  and  Z 
slowly  propagated  to  the  bottom. 

Finally,  two  sets  of  parameters  were  found  which  produced  well-defined  annual  cycles. 
One  without  a  spring  bloom  =  0.5)  was  characteristic  of  the  annual  cycle  in  the  subarctic 
North  Pacific  Ocean  in  the  vicinity  of  the  former  ocean  station  P.  (see  Fig.  7).  The 
other  set  of  parameters,  with  a  lower  initial  slope  of  the  P-I  fight  response  curve  (/?  = 
1.5),  produced  a  short  spring  bloom,  apparently  because  of  a  much  smaller  overwintering 
zooplankton  population,  less  able  to  crop  the  increase  in  phytoplankton  growth  as  the 
mixed  layer  shoals,  (see  Fig.  8).  (In  both  figures  the  initial  oscillations  result  from  the 

initial  conditions.) 

The  NPZ  model  coupled  with  the  MY  mixed  layer  model  (Figs.  9,  10)  reveals  that 
the  biological  variables  behave  almost  exactly  like  in  the  CONV  models,  given  that  the 
annual  cycle  of  the  mixed  layers  in  the  two  models  is  different.  In  MY  the  mixed  layer 
reached  the  bottom  of  the  model  at  80  m  each  winter  and  the  summer  mixed  layer  was 
not  as  shallow. 

The  original  form  of  CONV  indeed  had  a  30m  mixed  layer  depth.  For  the  particular 
forcing  function,  however,  the  PWP  and  MY  models  seemed  to  have  a  summer  MLD  of 
about  10m,  so  the  parameters  for  CONV  were  changed  to  use  this  value. 

An  NPZ  model  from  Peter  Franks,  adapted  to  include  detrital  material  (see  below), 
was  used  during  the  workshop  to  test  the  sensitivity  of  simple  models  to  detrital  sinking. 
The  remineralization  of  the  detrital  particles  is  treated  as  a  simple  first  order  process.  The 
NPZD  models  solutions,  shown  in  Figs.  11  and  12,  were  run  both  with  and  without  detrital 
sinking  (2  m/day),  and  the  total  nitrogen  content  of  the  surface  layer  is  2  mmol  m”®,  similar 
to  that  for  Bermuda.  Both  models  predict  a  strong  spring  bloom  in  both  phytoplankton 
and  zooplankton  with  very  low  nutrients  during  the  summer.  Adding  detrital  sinking 
modifies  the  model  solution  in  several  fashions:  surface  zooplankton  and  phytoplankton 
concentrations  are  lower  in  the  summer,  a  fall  bloom  occurs  in  the  sinking  case  when 
fall  mixing  entrains  nutrients  into  the  mixed  layer,  and  a  subsurface  nutrient  maximum 
develops  in  the  seasonal  thermocfine  due  to  detrital  remineralization.  The  subsurface 
nutrient  maximum  appears  to  be  an  unrealistic  feature  of  the  model,  resulting  from  too 
shallow  of  a  remineralization  scale. 

3.1.7  A  coupled  PWP  +  NPZD  model 

Hood,  Olson,  Flierl 

A  detritus  compartment  (D)  was  added  to  the  Flierl  and  Davis  (1993)  NPZ  biological 
model  (see  Table  III)  and  then  forced  with  a  modified  version  of  the  Price,  Weller  and 
Pinkel  (1986)  mixed  layer  model  (see  description  of  the  Hood/ Olson  PWP-NPZD  model 
in  Appendix  III.C.  for  details).  The  additional  detritus  compartment  can  be  considered  a 
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11.  NPZD/K-profile  model  solution  for  phytoplankton  (mM  N  m~^)  for  a  Sar¬ 
gasso  Sea  location  (64° W,  32° N)  with  no  sinking. 
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Fig.  12.  NPZD/K-profile  model  solution  for  phytoplankton  (mM  N  m“®)  for  a  Sar¬ 
gasso  Sea  location  (64° W,  32°N)  with  sinking  (2  m  d“^). 


Table  III.  Hood/Olson  NPZD  equations,  variables,  and  parameters. 


The  model: 

P  =  uptake  —  grazing  —  pd  x  P 
Z  =  ga  X  grazing  —  zd  x  Z 

N  =  (ae  —  ga)  x  grazing  —  uptake  +  Ae  xzdxZ  +  exD 
X)  =  (1  —  ae)  X  grazing  +  (1  —  Ae)  xzdxZ  +  pdxP  —  exD 

grazing  =  Rm  xPxZxXx{l  —  e  ^^) 
uptake  =  Vm  x  P  x  x  (j^) 


Variables  and  Coefficients: 


Description 

Symbol 

Value/Units 

Dissolved  nitrogen 

N 

pM 

Phytoplankton  nitrogen 

P 

pM 

Zooplankton  nitrogen 

Z 

pM 

Detritus  nitrogen 

D 

pM 

Light  intensity 

I 

dimensionless 

Max.  phytoplankton  growth  rate 

Vm 

0.3/day 

Max.  zooplankton  Grazing  rate 

Rm 

0.71  day 

Zooplankton  Death  rate 

zd 

0.0b ! day 

Detritus  remineralization  rate 

e 

0.2fday 

Phytoplankton  scenscence  rate 

pd 

O.Objday 

Growth  efficiency  for  Z 

ga 

0.35 

Assimilation  efficiency  for  Z 

ae 

0.70 

Assimilation  efficiency  for  higher 

Ae 

0.70 

trophic  levelsf 

Half  saturation  constant  for  N 

Ks 

O.lpM 

Half  saturation  constant  for  I 

P 

0.01 

Saturation  constant  for  grazing 

A 

3.5/P 

t  Growth  efficiency  for  higher  trophic  levels  set  at  zero 
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simple  “microbial  loop”,  where  the  detrital  remineralization  rate  (e)  is  the  rate  at  which 
bacteria  convert  detrital  Nitrogen  to  dissolved  Nitrogen. 

Before  running  the  coupled  model  we  (Steele,  Davis,  Fherl,  Olson  and  Hood)  chose  a 
set  of  coefficients  for  the  biological  model  (Table  III).  Most  of  these  are  accepted  values  for 
diatom/copepod  dominated  food  webs,  based  upon  observational  and  experimental  data. 
No  clear  consensus  was  reached  as  to  the  proper  magnitudes  of  the  zooplankton  death  rate 
(zd),  the  detritus  remineralization  rate  (e),  and  the  phytoplankton  senescence  rate  {pd). 

We  forced  the  PWP  model  with  idealized  functions  for  the  Atlantic  ocean  at  20°  N 
(Bleck  et  al.,  1989),  which  resulted  in  seasonal  variations  in  mixed  layer  depth  between 
10  and  100m.  Under  these  forcing  conditions  the  ecosystem  model  would  not  run  with 
the  coefficients  set  as  above,  primarily  because  the  phytoplankton  could  not  survive  the 
winter-time  deep  mixing  period.  In  order  to  overcome  this  problem  it  was  necessary  to 
increase  the  maximum  phytoplankton  growth  rate  (Vmax'j^  and  to  lower  the  half  saturation 
constant  for  hght,  (/3). 

The  coupled  model  is  sensitive  to  the  grazing  formulation.  With  simple  Ivlev  grazing 
we  were  unable  to  find  a  set  of  model  coefficients  that  gave  a  stable  ecosystem  over  the 
entire  range  of  physical  variability.  We  therefore  used  a  modified  version  of  the  Ivlev  equa¬ 
tion  (see  equations)  that  includes  an  additional  linear  dependence  on  P.  W^ith  this  grazing 
formulation  temporal  changes  in  phytoplankton  and  zooplankton  biomass  are  much  more 
tightly  coupled.  This  damps  oscillations  and  helps  prevent  runaway  growth  of  phytoplank¬ 
ton  in  the  summer  and  extinction  in  the  winter.  There  are  other,  possibly  better  ways 
to  prevent  runaway  growth  and  extinction.  Some  participants  in  the  workshop  recom¬ 
mended  using  the  Ivlev  formulation  with  a  minimum  grazing  threshold  (Moisan).  Steele 
recommended  using  a  formulation  of  the  form  grazing  =  i2Tn(p2/(l  +P^)),  which  results 
in  reduced  grazing  pressure  at  very  low  phytoplankton  populations.  The  incorporation 
of  an  additional  phytoplankton  species  with  different  growth  characteristics  might  be  the 
most  ecologically  realistic  way  to  solve  the  problem  (e.g.  include  one  species  adapted  to 
summer-time  conditions  and  another  adapted  winter-time  conditions). 

Previous  work  by  Hood  and  Olson  has  revealed  that  the  diffusion  subroutine  m  PWP 
is  non  conservative.  We  therefore  recommend  incorporation  of  a  better  numerical  scheme 
(e.g.  Crank-Nicholson). 

3.2  Food  Webs 

Armstrong  (Chair),  BoUens,  Frost,  Landry,  Landsteiner,  Moisan 
3.2.1  Motivation 

Our  overall  objective  was  to  construct  a  general  model  of  a  pelagic  food  web  and  to 
embed  this  model  in  an  appropriate  1-dimensional  physical  model.  Since  several  members 
of  the  group  have  extensive  experience  modelling  the  subarctic  Pacific  (e.g..  Frost  1987; 
Miller  et  al.  1991),  we  decided  to  develop  this  model  with  expHcit  reference  to  OWS  P. 
Two  characteristic  biogeochemical  features  of  this  region  are  that  nitrate  remains  high 
(>5-7  mmol  m~*)  and  chlorophyll  remains  low  (<1  mg  chi  a  m“®)  throughout  the  year 
(Miller  et  al.  1991).  In  addition,  the  phytoplankton  are  dominated  by  small  (<20  pm) 
cells  (Booth  1988),  but  blooms  of  large  cells  can  be  caused  by  addition  of  iron  to  natural 
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cultures  (Martin  et  al.  1991).  We  wished  to  see  under  what  circumstances  a  general  pelagic 
food  web  model  would  reproduce  these  features. 

3.2.2  Model  formulation 

The  biological  model  was  based  on  that  of  Frost  (1987);  that  model  was  modified 
to  include  multiple  size  classes  of  both  predators  (zooplankton)  and  prey  (phytoplankton 
and  smafier  zooplankton).  The  model  could  potentially  support  four  size  classes  each  of 
phytoplankton  and  zooplankton. 

We  considered  four  types  of  food  webs  (Fig.  13).  The  particular  scenario  we  tested 
was  whether,  given  realistic  physical  forcings  for  Station  P,  a  multichain  food  web  model 
would  collapse  to  one  containing  only  small  phytoplankton  species  while  allowing  nitrate 
to  remain  high,  or  whether  all  four  chains  would  persist,  producing  lower  mtrate  values 
than  are  typically  observed.  If  the  latter  happened,  we  would  then  test  whether  adding  a 
parameterized  form  of  iron  Hmitation  (Armstrong,  in  press)  would  make  the  model  behave 

more  hke  the  subarctic  Pacific.  .  ,  , 

Nominal  size  classes  for  phytoplankton  were  set  as  1,  4,  16,  and  64  fj,m  eqmvalent 
spherical  diameter;  zooplankton  size  classes  were  offset  from  those  of  phytoplankton  by  one 
length  unit  at  4,  16,  64,  and  256  fim.  Rate  constants  for  phytoplankton  and  zooplankton 
in  different  size  classes  were  assumed  to  follow  allometric  (power  law)  relationships  so  that 
the  rate  constant  Ri  for  the  ith  size  class  is  related  to  that  of  the  next  smaller  size  class 

ft  =  ft-i  (  Li  /  ij_.  (1) 


where  Li  is  the  characteristic  size  (length)  of  size  class  i  and  is  the  aUometric  con¬ 
stant  associated  with  rate  process  R  (Moloney  and  Field  1989,  1991).  Mesozooplankton 
(predominantly  copepods  in  the  genus  Neocalanus  were  included  as  upper  level  predators, 
consuming  all  other  elements  of  the  food  web  except  the  smallest  algal  size  class.  Ni¬ 
trogenous  substrate  (“nitrate”)  was  assumed  to  be  conserved  in  the  system,  with  instant 
remineralization  to  nitrate  upon  death  of  one  of  the  hving  components.  Bacteria  were 
excluded  from  consideration  since  their  inclusion  would  have  required  differentiating  am¬ 
monium  from  nitrate  so  that  the  balance  of  nutrient  regeneration  (as  ammonium),  export 
(as  one  or  more  forms  of  organic  nitrogen),  and  import  (as  nitrate  from  depth)  could 
be  expHcitly  modelled.  The  defining  equations  for  the  growth  of  phytoplankton  {Pi)  and 

zooplankton  {Zj)  are 


dPi 

dt 


—  Pi 


E  Z,^^Hj(F,)-XiPi 


Zj  eating  Pi 


and 


dZj 

dt 


Zj  jiHj(Fj)  -  E  Zk^SkiFt)  -  A,  Zj 


Zh  eating  Zj 


(2) 


(3) 


In  these  equations  Fj  is  the  food  available  to  zooplankton  size  class  j;  it  is  calculated 
as  the  unweighted  sum  of  the  densities  (measured  as  nitrogen  equivalents)  of  all  size  classes 
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of  phytoplankton  and  zooplankton  preyed  upon  by  zooplankton  size  class  j;  this  set  differs 
among  food  webs  (see  Fig.  13).  The  total  harvest  rate  Hj  of  food  by  zooplankton  size 
class  j  was  taken  to  be 


H, 


Fj  —  Fthre»h,j 


max,j 


Kztj  "i"  Fj  —  FthreahJ 


(4) 


where  Hmaxj  is  the  maximum  harvest  rate  [mmol  prey  (mmol  Z)-^  KzJ  is  the 

half-saturation  constant  for  zooplankton  feeding,  and  Fthreshj  is  a  threshold  food  level  at 
which  feeding  stops.  In  these  equations,  phytoplankton  growth  rate  for  size  class  i  {fn)  is 

also  given  as 


l^'TnaXfi 


N 

KN,i  +  N' 


(5) 


where  N  is  nitrate  concentration  (mmol  m"^),  Hmax,i  is  maximum  growth  rate  (d”^),  /(/) 
is  the  fraction  of  maximum  growth  allowed  at  Hght  level  I  (see  pp.  52-54  of  Frost  1987 
for  details),  and  KN,i  is  the  half-saturation  constant  for  nitrate  uptake.  Seasonal  data  on 
mesozooplankton  abundance  were  combined  with  data  on  Neocalanus  feeding  rates  (Frost 
1987)  to  supply  loss  rates  Aj  and  Xj  according  to  Equations  23  and  24  of  Frost  (1987). 

Total  nitrogen  {N  -f-  SPi  -f-  T.Zj)  was  taken  to  be  11.5  mmol  m^  Simulation  results 
were  generated  using  physical  forcings  appropriate  for  Station  P  (Fig.  2  of  Frost  1987). 
Values  for  critical  parameters  are  Usted  in  Table  IV;  all  concentration  units  are  nitrogen 

equivalents. 


Table  IV.  Model  parameters  used  in  the  food  web  simulations. 


Parameter 

Value  for  smallest  size  class 

AUometric  coefficient 

f^maXji 

1.4  d-^  @  12°C 

0,  -1 

Hmax,j 

2.0  d-i  @  5°C 

-0.25 

KN,i 

0.2  mmol  m^O.l 

0.5 

Kz,j 

0.13  mmol  m^  (10  fig  Chi  a  1“^) 

0.5 

Fthresh,j 

0.075  mmol  m®  (5  fig  C  1“^) 

0 

7i 

0.4 

0 

3.2.3  Results 

Model  runs  were  made  with  four  different  trophic  structures  (Fig.  13).  Initial  condi¬ 
tions  were  equal  biomass  concentrations  in  each  size  class  of  phytoplankton  and  each  size 
of  consumer.  Each  version  of  the  model  was  run  for  several  years  until  a  steady  annual 
cycle  was  obtained.  The  results  illustrated  below  are  steady  annual  cycles. 

Model  runs  with  trophic  structures  1-3  predicted  blooms  of  large-sized  phytoplankton. 
Larger  sizes  of  phytoplankton  bloomed  because  the  disparity  between  specific  growth  rates 


28 


of  phytoplankton  and  those  of  their  chief  consumers  increased  with  size.  The  timing  of 
phytoplankton  blooms  depended  on  trophic  structure.  With  trophic  structure  1  a  rektively 
small,  brief  bloom  of  large  cells  (P4)  occurred  at  the  time  of  spring  restratification  (Fig.  14); 
after  the  bloom  phytoplankton  stock  remained  low  and  nutrient  concentration  was  high. 
The  model  with  trophic  structure  2  produced  a  succession  of  blooms  -  first  a  brief  bloom 
of  P2  then  a  brief  bloom  of  PI,  and  finally  a  bloom  of  P3  that  persisted  throughout  the 
summer  and  early  fall  (Fig.  15).  In  this  case  N  was  depleted  in  summer.  Trophic  structure 
3  produced  an  annual  cycle  similar  to  that  for  trophic  structure  1,  namely,  a  smaU  bloom 
at  the  time  of  spring  restratification  followed  by  a  low  chlorophyll-high  N  condition  for  the 

rest  of  the  year  (Fig.  16). 

In  runs  of  the  model  with  trophic  structures  1-3  the  magnitude  of  phytoplankton 
blooms  was  increased  by  any  process  that  reduced  the  growth  rates  of  consumers  relative  to 
those  of  their  prey  organisms  (e.g.,  an  additional  mortality  rate  imposed  on  all  consumers, 
either  as  by  a  constant  daily  mortality,  or  predation  by  mesozooplankton).  The  magnitude 
of  phytoplankton  blooms  was  reduced  somewhat  by  decreasing  the  allometric  constant 
for  half-saturation  of  ingestion  by  consumers  (e.g.,  ^  for  Kz,i  decreased  to  0.25)  or  by 
introducing  weak  size-dependence  of  phytoplankton  maximum  specific  growth  rate  (e.g.,  ^ 
=  -0.25  for  /Xmax  i)-  However,  blooms  of  large  phytoplankton  cells  could  be  fully  suppressed 
only  under  strong  size-  dependence  of  phytoplankton  maximum  specific  growth  rate,  as 
might  pertain  under  iron  limitation.  For  example,  with  /3  =  -1.0  for  the  versions  of 

the  model  with  trophic  structures  1-3  produced  the  year-round  condition  of  low  chlorophyU- 
high  N  observed  in  the  open  subarctic  Pacific.  The  result  for  trophic  structure  1  is  shown 
in  Fig.  17.  In  these  runs  of  the  model,  one  or  two  of  the  larger  size  classes  of  phytoplankton 
(P3,  P4)  were  greatly  reduced  in  abundance  if  not  entirely  lost. 

The  model  results  suggest  that  strong  size-dependent  limitation  of  phytoplankton  spe¬ 
cific  growth  rate  is  required  to  reproduce  the  low  chlorophyU-lugh  N  condition  observed 
in  the  open  subarctic  Pacific.  However,  the  model  with  trophic  structure  4  reproduced 
the  condition  without  imposed  Hmitation  on  growth  rate  of  the  large  size  classes  of  phyto¬ 
plankton  (Fig.  18).  In  the  steady  annual  cycle  the  trophic  structure  reduced  to  a  simple 
Hnear  food  chain  consisting  of  PI  and  two  size  classes  of  consumers  (Z2  and  Z4).  Very 
similar  results  were  also  obtained  when  reasonable  levels  of  mortality  (either  a  constant 
daily  mortality  or  predation  by  mesozooplankton)  were  imposed  on  all  consumers.  How¬ 
ever,  even  with  trophic  structure  4  blooms  of  very  large  phytoplankters  might  be  predicted 
if  the  model  was  extended  to  include  one  or  two  larger  size  classes  of  phytoplankton  and 
consumers.  Although  this  trophic  structure  seems  most  realistic  (i.e.,  a  consumer  can  in¬ 
gest  any  trophic  type  of  prey  -  phytoplankton  or  zooplankton  -  in  the  size  classes  utilizable 
by  it),  an  elaboration  of  this  trophic  structure  might  be  to  assign  to  consumers  variable 
ingestion  efficiencies  on  different  sizes  of  prey. 

Predictions  of  the  model  with  trophic  structures  1-3  were  not  sensitive  to  the  level 
of  phytoplankton  production.  For  example,  increasing  a  (the  imtial  slope  of  the  photo¬ 
synthesis  vs.  irradiance  response  of  the  phytoplankton)  from  0.75  to  5.0  approximately 
doubled  the  annual  phytoplankton  production,  but  the  seasonal  cycles  of  abundance  of 
phytoplankton  and  consumers  were  essentially  the  same  as  described  above.  In  contrast, 
the  same  increase  in  a  in  the  model  with  trophic  structure  4  led  to  a  large  summer  bloom 
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of  PI  (not  a  feature  observed  in  the  open  subarctic  Pacific),  with  attendant  depletion  of  N, 
but  large  sizes  of  phytoplankton  still  did  not  persist.  The  low  chlorophyU-high  N  condition 

was  recovered  only  by  imposing  heavy  mortality  on  Z4. 

Perhaps  the  most  surprising  outcome  of  these  modelling  studies  was  that  introduction 
of  tropHc  structure  did  not  result  inevitably  in  predicted  annual  cycles  with  blooms  of 
large  phytoplankton  size  classes.  Indeed,  predictions  included  ’’balanced”  annual  cycles 
and  annual  cycles  with  blooms  of  the  smallest  size  class  of  phytoplankton.  These  few  ex¬ 
amples  are  sufficient  to  demonstrate  that  predictions  of  a  food  web  model  can  be  strongly 
influenced  not  only  by  the  specific  nature  of  the  size-dependence  of  physiological  processes 
modelled,  but  also  by  the  trophic  structure  assumed,  and  the  overall  level  of  biological 
producti^fity.  Knowledge  of  trophic  interactions  must  be  greatly  improved  before  a  partic¬ 
ular  trophic  structure  can  be  advocated  in  pelagic  foodweb  models  for  the  open  subarctic 
Pacific  or  other  pelagic  ecosystems. 


3.3  Structured  Populations 

Davis,  Olson  (Chair),  Pascual,  Steele 
3.3.1  Introduction 

The  idea  of  this  working  group  is  to  explore  means  of  taking  structure  in  a  population 
into  account  in  coupled  biological/physical  models.  Here  structure  can  be  defined  as  vari¬ 
ations  in  individuals  within  a  population  which  become  expressed  in  terms  of  population 
density,  n,  changes  in  either  space  and/or  time.  The  most  general  form  can  be  expressed 
in  terms  of  a  conservation  statement  for  the  population  which  is  expressed  in  an  expansion 
of  the  total  variation  of  n, 


^  dn  ,  dn  j  dn  , 
dn(t,s,w,£)  =  ^  ^ 


where  the  variables  in  the  expansion  are  time  (t),  hfe-stage  (i),  body-weight  (»),  and  space 
(x).  The  p  term  expresses  mortality  in  the  population.  The  total  change  in  time  of  the 

population  is  then  given  by 


dn 

dt 


dn  _  dn-  ,  r,-idn  ^  .  dn  _ 
-rr-  +  V  • +D  —  +5'^  — 
dt  dx  ds  dw 


-fin. 


Fluid  dynamicists  will  recognize  the  material  derivative  in  the  first  two  partial  expres¬ 
sions  although  here  one  might  also  include  mobile  behavior  on  the  part  of  the  population. 
The  age  or  stage  term  gives  the  equation  the  form  of  McKendrick-von  Foerster  (McK- 
endrick,  1926;  von  Foerster,  1959).  In  their  case  D  =  1  and  the  dependence  is  on  age. 
If  the  population  has  structure  based  on  life  stages  which  are  not  linear  in  time  then  D 
becomes  the  stage  duration  function  as  discussed  by  Corkett  and  McLaren  (1970)  and 
McLaren  (1978).  The  metabolic  term  is  most  commonly  tied  to  the  size  or  weight  of  the 
individual  (Sinko  and  Streifer,  1967;  de  Roos  et  a/.,  1992).  These  choices  of  expansion  vari¬ 
ables  are  not  unique.  For  example,  the  dependence  of  a  population  on  its  genetic  structure 
might  be  considered  (Webb,  1981).  Finally,  it  is  worth  distinguishing  the  changes  in  the 
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properties  of  the  equation  if  fi  =  K^),  the  equation  is  nonhnear.  Of  coarse,  it  is 
assumed  that  in  most  cases  fi  wiU  depend  on  the  expansion  variables,  c.f.  fx  - 

It  is  worthwhile  commenting  on  the  basic  problem  poised  by  the  McKendrick-von 
Foerster  equations  and  their  method  of  solution.  If  for  the  moment  interest  is  restricted 
to  the  case  with  D  =  1  and  ignoring  the  spatial  terms  the  problem  becomes 


dn/dt  +  dn/da  =  -fin. 

This  problem  has  both  initial  conditions,  n(a,0)  =  and  a  boundary  condition  based 
on  reproduction  in  the  population.  The  latter  can  be  written  in  terms  of  an  integral  of  a 
fecundity  function  (F)  over  age. 


F(a,w,jF)  do- 


Here  fecundity  is  assumed  to  be  a  function  of  age,  physiological  variables  (u;),  and  an 
effective  sex  ratio  for  the  population  (7f)-  Notice  that  the  initial  condition  and  the  birth 
rate  at  t  =  0  do  not  necessarily  match.  In  this  case  the  solution  will  have  initial  transients 

which  must  be  solved  for. 

There  are  various  means  of  seeking  solutions  to  the  McKendrick-von  Foerster  equation 
(c.f.  Webb,  1985,  Murray,  1989).  For  now  consider  the  hnear  problem,  fi  /  /^(^)-  One 
of  the  most  straight  forward  at  least  in  concept  is  to  introduce  a  cohort  function  for  a 
population  born  at  a  given  time,  i.e.  nc{t)  =  n{ac  +  t)  for  t  >  tc-  The  trajectory  of  the 
age  of  a  cohort  in  time  is  a  straight  hne  as  long  as  D  =  1  and  forms  a  characteristic  for  the 
solution  for  any  D.  The  problem  can  then  be  solved  for  a  individual  cohort  by  integration 
along  the  characteristic  as  shown  in  Fig.  19.  This  solution  for  a  given  cohort  is  then 


nc{t) 


nc{tc)exp[- 


df. 


This  can  then  be  extended  to  an  entire  range  of  cohorts  to  given  the  general  solution 


n{a,t) 


'  n(0,t  -  Oc)ea;p[-  /„“  t  >  a 

n{a  -  t,Q)exp[-  p{t')di%t  <  a 


For  simple  choices  of  p  these  integrals  can  be  solved.  Another  method  for  the  long 
term  hmit  on  the  population,  t  »  a^ax,  is  to  make  use  of  similarity  solutions  of  the 

formn(a,i)  =  n{a)e^*.  Upon  substitution  back  into  the  original  equation  and  some  work 

on  the  birth  condition  this  leads  to  viable  solutions.  Examples  of  similarity  solutions  are 
given  by  Murray  (1989). 

Our  goal  is  to  explore  structured  models  in  terms  of  some  concrete  examples  for  manne 
organisms  and  what  they  add  to  the  dynamics  of  these  populations.  In  particular  the  case 
of  copepods  with  their  multi-stage  structured  hfe-history  wiU  be  considered.  This  report 
presents  the  basic  copepod  model  formulated  at  the  workshop.  Then  this  basic  model 
is  embedded  into  an  NPZ  system.  We  briefly  present  directions  to  be  pursued  after  the 

workshop. 
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Fig.  19.  Dependence  of  stage,  5,  upon  time,  t,  for  an  organism  born  at  time  to.  For  a 
stage  duration,  D  =  1,  the  organism  is  described  by  a  linear  age,  a,  function.  For 
organisms  such  as  crustaceans,  which  follow  a  pattern  of  descrete  stages  [D  ^  1), 
the  stage  curve  will  in  general  take  some  other  form  and  fall  to  the  left  of  the 
hnear  age  relationship.  See  the  text  for  further  discusion. 
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3.3.2  Basic  copepod  model 

We  start  by  discretizing  the  cohorts  such  that  there  will  be  a  separate  equation  for  each 
cohort  as  identified  by  a  subscript  i.  The  model  then,  is  based  on  three  equations  governing 
the  weight  of  the  average  member  of  a  cohort,  the  number  of  individuals  m  a  cohort,  and  a 
diagnostic  relationship  for  age/stage  to  time  given  a  cohort’s  imtial  condition.  For  a  general 
preLtation  of  this  approach  see  de  Roos  et  al.  (1992).  The  basic  assumption  made  is  a 
Lphng  between  the  mean  condition  of  the  cohort  as  governed  by  the  mean  weight  of 
the  animals  in  a  cohort  of  size  m  through  the  influence  of  a  consumption  function, 
where  P  is  a  resource  and  T,  the  environmental  temperature,  and  mortahty, 

General  forms  for  these  can  be  written  as 

wi  =  fi{P,T)wi 


and 


lii  =  -mi{wi,T)niWi. 

The  basic  problem  then  is  to  find  appropriate  forms  for  mi{wi,T)  and 

relationship  for  the  transition  of  the  population  through  the  stages  which  make  up  the  hfe 

history  of  the  organism. 

As  an  example,  consider  a  popnlation  of  copepods  which  have  grorrth  curves  as  shown 
in  Fig.  20.  The  ultimate  growth  rate  is  set  by  the  exponential  growth  function  of  Huntley 
and  Lopez  (1992)  for  copepods.  An  individual  species  adapted  to  some  region  is  assume 
to  have  a  more  restricted  growth  function  centered  on  some  optimal  temperature,  2o.  Ihe 
weight  equation  is  then  of  the  form 

Wi  =  R„e°‘-^exp{-^{S){T-Tof}il-e-^^)wi 

where  the  (1  -  e-^^)  term  is  the  Ivlev  parameterization  for  resource  hmitation  and  S  is 
the  stage  of  the  copepod.  The  equation  for  the  number  of  individuals  in  a  cohort  is  given 

by 


rii  =  —[mo  +  Tni{wo{S,t)  —  Wi{S,t))/wo{S,t)]ni 

where  Wo  is  the  optimal  weight  at  stage  for  the  species  of  interest.  Tins  optimal  weight  can 
be  found  by  assuming  food  is  not  hmiting,  P  oo,  and  T  =  To.  Under  these  assumptions 


Wo  =  Woexp{Roe°‘'^'^t}. 

Finally,  we  need  a  relation  for  the  stage  of  the  copepod  as  a  function  of  time.  Here 
we  take  the  stage  duration  (D)  equation  of  form 


D  =  ai{T  —  as)^ 

after  McLaren  (1978).  The  stage  of  the  copepod  is  then  governed  by  the  modulus  of 


S  =  S  +dt/D{t) 
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Fig.  20.  Copepod  growth  as  a  function  of  environmental  temperature,  the  soHd  curve 
is  the  exponential  relationship  found  by  Huntley  and  Lopez  (1992)  using  a  broad 
compilation  of  growth  data.  The  dashed  curves  show  the  adaptive  modification 
of  this  curve  for  copepods  with  an  optimal  temperature  of  To  and  environmental 
range  of  half  width  1/ \/^.  A  hypothetical  limitation  on  growth  rate  is  shown  as 
a  horizon  below  which  the  population  becomes  nonviable. 
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where  S  has  a  value  of  one  for  eggs,  2-7  for  naupliar  stages  I- VI,  and  8-13  for  the  copepodid 
stages,  the  last  stage,  13,  being  the  adult  copepod.  We  are  applying  these  equations 
to  the  copepod  Calanus  finmarchicus  and  have  also  incorporated  temperature-dependent 
diapause  (a  resting  stage)  in  the  fifth-stage  copepodid. 

3.3.3  Coupling  the  model  to  the  rest  of  the  ecosystem 

The  basic  model  in  this  work  pictures  a  system  consisting  of  nutrients,  TV,  phyto¬ 
plankton,  P,  and  zooplankton,  Z,  where  the  latter  is  assumed  to  consist  of  a  structured 
population.  All  of  these  variables  have  units  of  [pmol/kg].  A  total  population  of  zoo¬ 
plankton  is  divided  into  cohorts  such  that,  Zt  -  with  Zi  =  riiWi  where,  as 

before,  Wi  denotes  the  mean  individual  weight  ^pmol f  kg jindividual],  and  n,,  the  number 
of  individuals  of  cohort  i. 

The  total  balance  in  the  system  is  then  governed  by 

N  =  -aNP  -f  ^mi{wi,T)niWi 

i 

P  =  aNP  - 

i 

Zi  =  fi{P,T)niWi  -  mi{wi,T)niWi 

where  the  functions  mi{wi,T)  and  MP,T)  denote  the  mortality  and  growth  of  the  mean 
individual  in  each  cohort. 

We  are  investigating  if  and  how  incorporating  structure  modifies  the  dynamics  of 
the  NPZ  system.  We  are  presently  considering  the  physical  environment  only  through 
temporal  forcing  (i.e.,  the  yearly  temperature  cycle  and  changes  in  light  and  nutrient 
supply  as  affected  by  changes  in  mixed  layer  depth).  We  plan  to  use  this  model  in  our 
ongoing  biological-physical  modehng  studies  which  examine  mixed-layer  dynamics  as  well 
as  two  and  three-dimensional  circulations. 
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oxygen 

Description  -  He  gets  new  production  rates  at  Papa  of  50-100  gC/m2y  (02  su- 
persat  method)  and  10-60  gC/m2y  (N03  consumption  method). 

Kim,  J-.W.  1976.  A  generalized  bxilk  model  of  the  oceanic  mixed  layer.  J.  Phys. 
Ocean.,  6,  686-695. 
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Keywords  —  upper  ocean,  mixed  layer,  ID,  bulk  model,  wind  stress  scabng,  com¬ 
mented,  and  DOE 

Description  -  Related  to  the  mixed  layer  embedded  model  in  OSU’s  GCM? 

Kitaigorodskii,  S.A.  1981.  On  the  theory  of  the  surface-stress  induced  entrainment 
at  a  buoyancy  interface  (toward  interpretation  of  KP  and  KPA  experiments). 
Tellus,  33,  89-101. 

Keywords  —  upper  ocean,  mixed  layer,  ID,  bulk  models,  model  comparison,  Ri 
formulation,  wind  stress  scaling,  physical  basis,  commented,  and  DOE 
Description  -  Attempts  to  reconcile  the  data  of  Kato  and  Philhps  (69)  with  shear¬ 
generated  turbulence  ideas.  See  also  deSzoeke  and  Rhines,  76. 

Klein,  P.,  Coantic,  M.  1981.  A  numerical  study  of  turbulent  processes  in  the  marine 
upper  layers.  J.  Phys.  Ocean.,  11,  849-863. 

Keywords  -  upper  ocean,  mixed  layer,  ID,  turbulence  closure,  commented,  and 
DOE 

Description  -  They  use  Mellor’s  L3  model.  They  find  that  deepening  occurs  due 
to  the  shear,  rather  than  by  the  wind  directly. 

Klein,  P.,  Coste,  B.  1984.  Effects  of  wind-stress  variabihty  on  nutrient  transport  into 
the  mixed  layer.  Deep  Sea  Res.,  31(1),  21-37. 

Keywords  -  upper  ocean,  mixed  layer,  ID,  turbulence  closure,  new  production, 
flux  across  thermocline,  and  commented 

Description  -  They  use  Mellor’s  L2  model,  and  find  that  mixing  occurs  due  to 
shear,  rather  then  directly  by  the  wind.  They  find  that  nuts  mix  up  in  pulses, 
rather  than  steadily. 

Kraus,  E.B.,  Turner,  J.S.  1967.  A  one- dimensional  model  of  the  seasonal  thermochne 
II.  The  general  theory  and  its  consequences.  Tellus,  19,  98-105. 

Keywords  —  upper  ocean,  mixed  layer,  ID,  bulk  model,  wind  stress  scaling,  com¬ 
mented,  and  DOE 

Description  -  The  classic.  Later  improvements  were  made  in  the  dissipation  term 
(many  authors). 

KroU,  J.  1988.  Instability  of  a  mixed  later  model  and  the  generation  of  near-inertial 
motion.  Part  I:  Constant  mixed  layer  depth.  J.  Phys.  Ocean.,  July,  963. 
Keywords  —  upper  ocean,  mixed  layer,  and  near  inertial  motions 

Martin,  P.J.  1985.  Simulation  of  the  mixed  layer  at  OWS  November  and  Papa  with 
several  models.  J.  Geophys.  Res.,  90(C1),  903-916. 

Keywords  -  upper  ocean,  mixed  layer,  model  comparison,  commented,  and  DOE 
Description  -  See  the  technical  report. 

Martin,  P.J.  1986.  Testing  and  comparison  of  several  mixed-layer  models.  Naval 
Ocean  Research  and  Development  Activity  Tech  Rep,  143. 

Keywords  -  upper  ocean,  mixed  layer,  model  comparison,  commented,  and  DOE 
Description  -  Gives  computer  costs,  memory  constraints,  and  how  well  the  models 
simulate  the  observed  data,  both  on  long  and  short  time  scales.  1)  Data  from 
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FLIP  in  PWP  not  weU  simulated  by  the  turbulence  closure  models  -  they  must 
not  mix  in  the  mixed  layer  by  the  complete  mechanism.  2)  KT  models  have  that 
disturbing  ’’jump”  at  the  base  of  the  mixed  layer.  3)  The  various  turbulence 
closure  models  only  dilfer  in  their  Ri-crit  for  mixing.  They  are  mostly  shear 
models,  due  to  low  turbulence  transport  in  the  models. 

Meehl,  G.A.,  Washington,  W.M.  1985.  Sea  surface  temperatures  computed  by  a 
simple  ocean  mixed  layer  coupled  to  an  atmospheric  GCM.  J.  Phys.  Ocean.,  15, 
92-104. 

Keywords  -  upper  ocean,  AGCM,  coupled  A-0,  SST,  NCAR,  and  with  slab  mixed 
layer 

Description  -  Atmos.  GCM  with  a  slab  ocean  mixed  layer.  They  find  that  the 
SST  is  fairly  weU  represented  most  places,  but  to  high  in  eq.  upwelling  areas  (no 
advection  in  ocean).  Mixed  layer  is  of  a  fixed  depth,  50m. 

Mellor,  G.L.  1973.  Analytic  prediction  of  the  properties  of  stratified  planetary  surface 
layers.  J.  Atm.  Sci.,  30,  1061-1069. 

Keywords  -  upper  ocean,  mixed  layer,  ID,  TKE  model,  commented,  and  DOE 
Description  -  One  of  the  originals 

MeUor,  G.L.,  Yamada,  T.  1974.  A  hierarchy  of  turbulence  closure  models  for  planetary 
boundary  layers.  J.  Atmos.  Sci.,  31,  1791-1806. 

Keywords  -  upper  ocean,  mixed  layer,  ID,  turbulence  closure,  commented,  and 
DOE 

Description  -  Adds  the  simplification  scheme  to  get  all  the  different  levels  of 
models. 

Mellor,  G.L.,  Durbin,  P.A.  1975.  The  structure  and  dynamics  of  the  ocean  surface 
mixed  layer.  J.  Phys.  Ocean.,  5,  718-728. 

Keywords  —  upper  ocean,  mixed  layer,  ID,  turbulence  closure,  commented,  and 
DOE 

Description  -  ’’Comparisons  of  predictions  and  data  are  favorable.” 

Mellor,  G.L.  1979.  Retrospective  on  oceanic  boundary  layer  modeling  and  second 
moment  closure,  (in  manuscript). 

Keywords  -  upper  ocean,  mixed  layer,  ID,  turbulence  closure,  commented,  and 
DOE 

Description  -  Another  summary  of  the  physics,  reflecting  recent  discussion  in  the 
field. 

Mellor,  G.L.,  Yamada,  T.  1982.  Development  of  turbulence  closure  model  for  geo¬ 
physical  fluid  problems.  Rev.  Geophys.  Space  Phys.,  20(4),  851-875. 

Keywords  -  upper  ocean,  mixed  layer,  ID,  turbulence  closure,  commented,  and 

DOE 

Description  -  And  away  we  go. 

Miyakoda,  K.,  Rosati,  A.  1984.  The  variation  of  sea  surface  temperature  in  1976  and 
1977  2.  The  simulation  with  mixed  layer  models.  J.  Geophys.  Res.,  89(C4), 
6533-6542. 
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Keywords  —  upper  ocean,  mixed,  layer,  ID,  model  comparison,  bulk  model,  wind 
stress  scaling,  turbulence  closure,  and  commented 

Description  -  They  spread  out  ID  mixed  layer  models  over  the  Pacific  to  figure  out 
why  the  temperatures  were  anomalous  during  1977.  No  transport  (1  by  N  model, 
like  Glover’s).  Found  that  both  types  under-mixed  the  water  between  surface  and 
deep,  and  concluded  that  Mellor  level  2.0  was  not  good  enough.  Concluded,  hke 
Rosatti  (8X)  that  monthly  mean  wind  data  is  not  good  enough. 

Mourn,  J.N.,  CaldweU,  D.R.,  Paulson,  C.A.  1989.  Mixing  in  the  equatorial  surface 
layer  and  thermocHne.  J.  Geophys.  Res.,  94  (C2),  2005-2023. 

Keywords  -  upper  ocean,  equatorial  ocean,  mixed  layer,  turbulence,  and  dissipa¬ 
tion 

Musgrave,  D.L.,  Chou,  J.,  Jenkins,  W.J.  1988.  Application  of  a  model  of  upper-ocean 
physics  for  studying  seasonal  cycles  of  oxygen.  J.  Geophys.  Res.,  93,  15679- 
15700. 

Keywords  -  upper  ocean,  mixed  layer,  ID,  new  production,  dissolved  oxygen, 
bulk  model,  Ri  formulation,  commented,  and  DOE 

Description  -  Gets  new  production  rates  of  3-4  mol  02/m2y,  which  is  about  32 
gC /m2y,  as  opposed  to  50  (from  Jenkins  and  Goldman,  85).  Paper  also  interesting 
in  that  he  had  to  put  in  a  diffusivity  to  get  the  annual  SST  cycle  right.  P.Martin 
didn’t  have  to  do  that. 

Niiler,  P.P.,  Kraus,  E.B.  1977.  One- dimensional  models  of  the  upper  ocean.  In:  Mod- 
elling  and  Prediction  of  the  Upper  Layers  of  the  Ocean,  Kraus  (Ed.),  Pergammon 
Press,  Oxford,  143-172. 

Keywords  -  upper  ocean,  mixed  layer,  ID,  bulk  model,  wind  stress  scaling,  com¬ 
mented,  and  DOE 
Description  -  The  original. 

Pacanowski,  R.C.,  Philander,  S.G.H.  1981.  Parameterization  of  vertical  mixing  in 
numerical  models  of  tropical  oceans.  J.  Phys.  Ocean.,  11,  1443-1451. 

Keywords  —  upper  ocean,  mixed  layer,  thermochne  ventilation,  OGCM,  Ri  for¬ 
mulation,  commented,  and  DOE 

Description  -  Outhne  a  Ri-dependent  eddy  mixing  formulation.  These  vary  con¬ 
tinuously  with  Ri,  rather  than  being  turned  on  at  a  critical  Ri,  as  in  PIVP.  Also, 
they  are  gradient  Ri’s,  not  bulk  Ri’s.  PWP  was  unable  to  make  his  model  work 
with  just  Ri(grad).  Perhaps  the  difference  is  the  lower  resolution  of  this  model 
(5  boxes  in  upper  100m). 

Phillips,  O.M.  1966.  The  dynamics  of  the  upper  ocean,  Cambridge  Press,  London. 
Keywords  -  upper  ocean,  fluid  mechanics,  turbulence,  mixed  layer,  and  DOE 

Pollard,  R.T.,  Rhines,  P.B.,  Thompson,  R.O.R.Y.  1973.  The  deepening  of  the  wind- 
mixed  layer.  Geophys.  Fluid  Dyn.,  3,  381-404. 

Keywords  -  upper  ocean,  mixed  layer,  ID,  bulk  model,  Ri  formulation,  physical 
basis,  commented,  and  DOE 
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Description  -  Elegant  paper,  laying  out  lots  of  reasons  why  the  bulk  Ri  formula¬ 
tion  is  more  correct  than  the  wind,  stress  scahng.  See  also  deSzoeke  and  Rhines, 
76. 

Pollard,  R.T.,  Thomas,  K.J.H.  1989.  Vertical  circulation  revealed  by  diurnal  heating 
of  the  upper  ocean  in  late  winter,  part  I:  observations.  J.  Phys.  Ocean.,  19,  269. 
Keywords  —  upper  ocean,  mixed  layer,  convective  mixing,  data,  and  commented 
Description  -  They  observe  convective  mixing  in  the  surface  ocean. 

Price,  J.F.,  Mooers,  C.N.K.,  Van,  Leer.  1978.  Observations  and  simulation  of  storm- 
induced  mixed-layer  deepening.  J.  Phys.  Ocean.,  8,  582-599. 

Keywords  -  upper  ocean,  mixed  layer,  ID,  bulk  model,  model  comparison,  Ri 
formulation,  wind  stress  scaling,  commented,  and  DOE 

Description  -  Argues  strongly  that  the  Ri  scaling  is  superior  to  the  wind  scaling, 
with  some  observations  and  comparisons  with  the  model.  Finds  that  the  empirical 
constants  needed  for  wind  scaling  change  under  different  conditions,  whereas  for 
the  Ri  scaling,  they  are  pretty  constant. 

Price,  J.F.,  Weller,  R.A.,  Pinkie,  R.  1986.  Diurnal  cycling:  observations  and  models 
of  the  upper  ocean  response  to  diurnal  heating,  cooling  and  wind  mixing.  J. 
Geophys.  Res.,  91(C7),  8411-8427. 

Keywords  -  upper  ocean,  mixed  layer,  ID,  bulk  model,  Ri  formulation,  com¬ 
mented,  and  DOE 

Description  -  The  new  classic  Ri  model,  with  the  innovation  of  the  smoothing  of 
the  base  of  the  mixed  layer  using  a  Ri-gradient  formulation,  which  slows  it  down 
but  gives  a  more  realistic  water  column  profile  (see  Martin’s  technical  report). 
Also  some  great  high-resolution  diurnal  mixed  layer  current  meter  and  CTD  data 
from  FLIP,  which  shows  higher  mixing  in  the  ’’mixed”  layer  than  is  predicted  by 
the  turbrilence  closure  models  (see  again  Martin). 

Rosati,  A.,  Miyakoda,  K.  1988.  A  general  circulation  model  for  upper  ocean  simula¬ 
tion.  J.  Phys.  Ocean.,  18,  1601-1626. 

Keywords  -  upper  ocean,  mixed  layer,  OGCM,  Princeton  -  GFDL,  turbulence 
closure,  commented,  and  DOE 

Description  -  OGCM  skewed  toward  upper  ocean  detail.  6  levels  in  upper  70m. 
Mellor  level  2.5  turbiilence  closure.  12-hr  average  winds  worked  better  than 
monthly  averages  (evap  heat  flux  underestimated  with  monthly  ave.).  Heat  fluxes 
calcrdated  using  biolk  formulae.  No  mention  of  sea  ice. 

Sarachik,  E.S.  1978.  Boundary  layers  on  both  sides  of  the  tropical  ocean  surface.  In: 
The  Global  Weather  Experiment:  FINE  workshop  proceedings,  Nova  /  N.Y.I.T. 
Univ.  Press. 

Keywords  —  upper  ocean,  mixed  layer,  and  SST 

Sarmiento,  J.L.,  Fasham,  M.J.R.,  Siegenthaler,  U.,  Najjar,  R.,  ToggweUer,  J.R.  1989. 
Models  of  chemical  cychng  in  the  oceans:  progress  report  II.  Ocean  Tracers  Lab¬ 
oratory  Technical  Report  No.  6. 
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Keywords  -  upper  ocean,  mixed  layer,  new  production,  surface  ocean  bio  model, 
OGCM,  DOC,  commented,  and  DOE 

Description  -  Uses  the  model  of  Toggweiler  (JGR  94  (C6),  8217-8265,  ’89).  Also 
describes  the  nutrient-restoring  model  used  to  argue  that  DOC  provides  an  in¬ 
visible  nitrogen  pool  that  is  advecting  around. 

Schopf,  P.S.,  Cane,  M.A.  1983.  On  equatorial  dynamics,  mixed  layer  physics,  and  sea 
surface  temperatures.  J.  Phys.  Ocean.,  13,  917-935. 

Keywords  -  upper  ocean,  equatorial  ocean,  mixed  layer,  El  Nino,  coupled  A-0, 
SST,  and  commented 

Description  -  Use  bulk  mixed  layer  physics  to  do  ENSO  stuff 

Spitzer,  W.S.,  Jenkins,  W.J.  1989.  Rates  of  vertical  mixing,  gas  exchange,  and  new 
production:  estimates  from  seasonal  gas  cycles  in  the  upper  ocean  near  Bermuda. 
J.  Mar.  Res.,  47,  169-196. 

Keywords  -  upper  ocean,  mixed  layer,  ID,  thermochne  ventilation,  new  produc¬ 
tion,  dissolved  oxygen,  gas  exchange,  commented,  thermochne,  diapycnal  mixing, 
and  DOE 

Description  -  He  gets  a  vertical  turbulent  diffusivity  of  about  lE-4  m2/s,  based 
on  the  temperature  history  and  the  argon  data.  Gets  new  production  rates  of 
about  4-6  mol  02/m2y  (that’s  about  50g  C/m2y,  consistent  with  Jenkins  and 
Goldman).  Also  estimates  rates  of  air  injection  and  of  gas  exchange. 

Stevenson,  J.W.  1979.  On  the  effect  of  dissipation  on  seasonal  thermochne  models. 
J.  Phys.  Ocean.,  9,  57-64. 

Keywords  —  upper  ocean,  mixed  layer,  ID,  bulk  model,  wind  stress  scahng,  model 
comparison,  dissipation,  commented,  and  DOE 

Description  -  Extensive  review  of  the  issue  of  dissipation  in  the  wind-driven  bulk 
mixed  layer  models.  Some  of  them  mix  aU  the  way  to  the  bottom,  eventually. 

Stevenson,  J.W.  1983.  The  seasonal  variation  of  the  surface  mixed-layer  response  to 
the  vertical  motions  of  hnear  rossby  waves.  J.  Phys.  Ocean.. 

Keywords  -  abstract  only,  upper  ocean,  mixed  layer,  and  Rossby  waves 

Tabata,  S.  1961.  Temporal  changes  of  sahnity,  temperature,  and  dissolved  oxygen 
content  of  the  water  at  station  ”P”  in  the  Northeast  Pacific  Ocean,  and  some  of 
their  determining  factors.  J.  Fish.  Res.  Bd.  Canada,  18(6),  1073-1124. 
Keywords  —  upper  ocean,  mixed  layer,  ID,  and  data 

Tabata,  S.,  Boston,  N.E.J.,  Boyce,  F.M.  1965.  The  relation  between  wind  speed 
and  summer  isothermal  surface  layer  of  water  ot  ocean  station  P  in  the  Eastern 
Subarctic  Pacific  Ocean.  J.  Geophys.  Res.,  70,  3867-3878. 

Keywords  —  upper  ocean,  mixed  layer,  ID,  and  data 

Therry,  G.,  Lacarrere,  P.  1983.  Improving  the  eddy  kinetic  energy  model  for  planetary 
boundary  layer  description.  Boundary -Layer  Meteorology,  25,  63-88. 

Keywords  -  upper  ocean,  mixed  layer,  ID,  turbulence  closure,  commented,  and 

DOE 
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Description  -  Compared  by  P.  Martin  to  other  types  of  models.  Very  similar 
to  MeUor’s  but  they  come  out  with  slightly  different  values  for  the  critical  Ri 
number,  on  which  the  response  of  the  models  depend. 

Thomas,  A.C.,  Strub,  P.T.  1989.  Large-scale  patterns  of  phytoplankton  pigment 
distribution  during  the  spring  transition  along  the  west  coast  of  North  America, 
(manuscript). 

Keywords  —  upper  ocean,  new  production,  wind  forcing,  and  mixed  layer 

Thomas,  K.J.H.  1989.  Vertical  circulation  revealed  by  diurnal  heating  of  the  upper 
ocean  in  late  winter,  part  II:  modehng.  J.  Phys.  Ocean. ^  19,  279. 

Keywords  -  upper  ocean,  convective  mixing,  and  mixed  layer,  ID 

Venrick,  E.L.,  McGowan,  J.A.,  Cayan,  D.R.,  Hayward,  T.L.  1987.  Chmate  and 
Chlorophyll  a:  long-term  trends  in  the  central  North  Pacific  Ocean.  Science, 
238,  70-72. 

Keywords  -  upper  ocean,  mixed  layer,  new  production,  and  wind  forcing 
Description  -  Reports  that  chlorophyll  concentrations  in  the  North  Pacific  have 
increased  over  the  last  20  years.  Correlates  that  with  in  increase  in  winter  winds 
and  a  decrease  in  SST. 

Woods,  J.D.  1980.  Diurnal  and  seasonal  variation  of  convection  in  the  wind-mixed 
layer  of  the  ocean.  Quarterly  J.  Royal  Met.  Soc.,  106,  379-394. 

Keywords  —  upper  ocean,  mixed  layer,  ID,  diurnal  cycle,  seasonal  cycle,  convective 
mixing,  and  commented 

Description  -  He  says  that  wind  is  only  important  for  mixing  during  the  day, 
and  convection  only  at  night.  Modeling  this  requires  a  diurnal  timestep,  but  he 
presents  a  parameterization  which  does  this  with  longer  time  stepping,  for  chmate 
models. 

Woods,  J.D.,  Barkmann,  W.  1986.  The  response  of  the  upper  ocean  to  solar  heating. 
I:  The  mixed  layer.  Quart.  J.  Roy.  Met.  Soc.,  112,  1-42. 

Keywords  —  upper  ocean,  mixed  layer,  ID,  optics,  commented,  and  DOE 
Description  -  1.5  deg  K  error  in  SST  results  from  inadequate  resolution  of  diurnal 
cycle  (relative  to  resolved  model,  not  relative  to  data).  Also,  1  deg  K  error  from 
-f-  1  Jerlov  scale  unit  in  turbidity,  about  the  uncertainty  (see  page  20). 

Woods,  J.D.,  Barkmann,  W.  1986.  A  lagrangial  mixed  layer  model  of  Atlantic  18C 
water  formation.  Nature,  319,  574-576. 

Keywords  —  upper  ocean,  mixed  layer,  ID,  and  lagrangian  model 

Zeman,  O.,  Lumley,  J.L.  1976.  Modehng  buoyancy  driven  mixed  layers.  J.  Atl.  Sci., 
33, 1974-1988. 

Keywords  -  upper  ocean,  mixed  layer,  ID,  turbulence  closure,  and  commented 
Description  -  Cited  by  P.  Martin  (tech  report)  as  being  criticism  of  MeUor  model’s 
transport  of  turbulence  (as  being  inadequate). 
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Appendix  III.  Descriptions  of  some  existing  1-D  coupled  biological-physical 
models 


A.  Archer  PAPA  model 

David  Archer,  University  of  Chicago 

A  one-dimensional  model  of  temperature,  salinity,  nutrients,  oxygen,  carbon,  and  ar¬ 
gon  chemistry  is  used  to  hindcast  the  annual  cycle  of  sea  surface  pC02  at  weathership  sta¬ 
tion  Papa  in  the  subarctic  Pacific  (50°N,  145°W),  based  on  recent  biological  and  chemical 
measurements  made  in  conjunction  with  the  SUPER  program.  Heat  fiuxes  are  calculated 
from  the  meteorological  time  series  data  from  the  Canadian  weathership  program.  The 
Price  et  al.  (1986)  model  is  used  for  predicting  mixed  layer  dynamics.  The  rate  of  new 
production  in  the  model  is  based  on  sediment  trap  data  (Welschmeyer,  unpubhshed)  and 
a  comparison  of  model  oxygen  and  argon  concentrations  with  the  recent  data  of  Emerson 
et  al.  (1991).  The  balances  of  nutrients  and  oxygen  in  the  permanent  halochne  require 
isopycnal  ventilation  on  a  time  frame  of  approximately  10  years;  this  estimate  is  consistent 
with  the  estimate  of  VanScoy  et  al.  (1991)  based  on  tritium  data  from  Geosecs  and  Long 
Lines  programs.  The  model  results  are  compared  with  the  5  year  time  series  data  presented 
by  Wong  et  al.  (1990).  The  model  appears  to  capture  the  mean  sea  surface  pC02  and 
the  magnitude  and  timing  of  the  seasonal  cycle.  The  data  however  contain  much  greater 
high  frequency  variation  than  the  model  results,  which  we  beheve  is  caused  by  patchiness 
in  the  horizontal  distribution  of  NOS.  The  model  pC02  sensitivity  to  the  chemistry  of  the 
upwelling  water  and  the  rate  of  biological  new  production  are  presented.  Although  the 
model  simulation  of  pC02  is  satisfactory,  this  study  illustrates  the  hmitations  of  modeling 
the  chemistry  of  the  upper  ocean  in  one  dimension.  The  slow  currents  and  horizontally 
homogeneous  sea  surface  in  the  subarctic  Pacific  make  Papa  one  of  the  best  available  can¬ 
didates  for  modehng  in  1-D.  In  spite  of  this,  a  1-D  formulation  is  inadequate  to  address 
the  chemistry  of  the  halocline  (a  crucial  lower  boundary  condition  to  the  mixed  layer)  and 
does  not  constrain  the  transport  of  the  nutrients  by  wind-driven  currents  in  the  mixed 
layer.  We  conclude  that  further  progress  in  modeling  the  upper  ocean  nutrient  and  carbon 
cycles  will  require  simulation  in  three  dimensions. 
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B.  PWPBIO_OPT;  A  Fully-Coupled  One-Dimensional  Physical-Bio-Optical 
Model 

John  Moisan,  Scripps  Institution  of  Oceanography 

1.  Model  Description 

The  model  is  a  one- dimensional  (vertical)  time-dependent  physical-bio-optical  model.  The 
physical  portion  includes  the  effects  of  vertical  velocity  (biological  and  advective),  back¬ 
ground  vertical  diffusion,  and  Price- Weller-Pinkel  mixed-layer  dynamics  (Price  et  al.,  1986; 
hereafter  PWP).  The  biological  portion  of  the  model  consists  of  a  general  iV-component 
food  web  which  can  be  easily  reconfigured  to  the  specifications  required  for  other  mixed- 
layer  bio-optical  modeling  programs.  As  an  example,  the  model  is  currently  configured  with 
jjine  component  food-web  that  includes  silicate,  nitrate,  ammonium,  two  phytoplankton 
size  fractions,  three  zooplankton  categories  and  a  detrital  pool  (Moisan,  1992).  Transfers 
within  this  food  web  occur  through  nutrient  uptake  by  phytoplankton,  differential  grazing 
by  zooplankton  and  nutrient  recychng.  Finally,  coupled  to  the  general  food  web  model  is  a 
modified  version  of  the  subsurface  spectral  irradiance  model  of  Sathyendranath  and  Platt 

(1988). 

2.  Physical  Model 

The  governing  equations  which  are  used  in  this  model  are  the  one-dimensional  primi¬ 
tive  equations: 


du  du  .  d 


dv  dv  d 


dT  dT  d  dT\ 


8S  8S  d  ( ^  as\ 
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bu - r— 

p  oz 


J~ y  bv  - 
p  oz 


pc  \  dz 


+  PHAR 


^  dE 


(1) 

(2) 

(3) 

(4) 


where  the  terms  and  parameters  in  equations  1-4  are  defined  in  Table  1. 

Equations  1  and  2  express  the  u,  v  momentum  balances,  respectively.  The  time 
evolution  of  the  temperature,  T{z,t),  and  salinity,  S{z,t),  fields,  which  govern  the  density, 
p{z,t),  are  given  by  equations  3  and  4,  respectively. 

The  vertical  velocity  at  the  sea  surface,  the  top  of  the  model  domain,  is  set  to  zero. 


w{z  =  0  m,t)  =  0. 

The  vertical  velocity  at  the  bottom  of  the  domain,  Zmax^  is  defined  as. 

w(z  =  •Zmo*j^)  — 
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Table  1.  Definition  of  Terms  in  the  Primitive  Equations  For  the  One- 
Dimensional  Model 


Terms 
u{z^  t) 
v{z,t) 
w{z,t) 
Tiz,t) 
S{z,t) 
p{z,t) 


f 

9 

8 

K,{z) 

I{z) 

PHAR{z) 

Eiz) 


Definition  _ _ 

horizontal  velocity  in  the  x-direction 

horizontal  velocity  in  the  y-direction 

vertical  velocity 

temperature 

salinity 

density 

wind  stress  in  the  x-direction 
wind  stress  in  the  y-direction 
specific  heat  of  seawater 
Coriolis  parameter 
acceleration  of  gravity 
linear  drag  coefficient 
background  diffusion  term 
climatological  forcing  terms 
subsurface  solar  irradiance 
photosyntheticaUy  absorbed  radiation 
evaporation  minus  precipitation 
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where  f{t)  represents  a  time-varying  or  constant  vertical  velocity.  The  interior  vertical 
velocities  are  then  obtained  by  assuming  a  hnear  dependence  in  the  vertical  between  the 
surface  and  the  bottom  of  the  model  domain, 

dw{z,t)  _  w(z  =  Zmaxji)  ~  —  0 

dz  Zmax 

The  interior  vertical  velocities  are  then  obtained  from: 


w{z,t)  = 


w(^Z  —  Zfnax^^') 


where  z  is  the  depth.  This  method  assumes  a  constant  horizontal  divergence  throughout 
the  model  domain.  Finally,  a  background  vertical  diffusion  term,  Kz{z),  is  included  in  the 
model  as  a  depth-dependent  term.  By  doing  this,  diffusion  can  be  increased  in  regions  of 
the  domain  where  wind  mixing  may  not  play  an  important  role  but  mixing  due  to  other 
processes  is  known  to  occur. 


3.  PWP  Mixed-Layer  Dynamics 

The  model  incorporates  PWP  mixed-layer  dynamics  (Price  et  al.,  1986).  Vertical 
mixing  occurs  in  the  model  in  order  to  meet  the  following  three  criteria  of  the  PWP 
mixing  scheme: 

1)  The  first  mixing  criteria  simulates  free  convective  mixing  in  the  water  column  in 
order  to  achieve  static  stability.  This  requires  that 


(5) 


2)  The  second  criteria  simulates  the  entrainment  of  the  mixed-layer.  Termed  mixed- 
layer  stability,  this  entrainment  process  is  accomplished  by  a  relaxation  of  the  bulk  Richard¬ 
son  number 


Rb 


gAph 

poiAvy 


>  0.65. 


(6) 


3)  Finally,  mixing  due  to  the  effect  of  shear  flow  instabihty  is  simulated  by  a  relaxation 
of  a  gradient  Richardson  number 


R 


9 


gdpidz 

p^idVldzf 


>  0.25, 


(7) 


The  forcing  of  the  model  is  accomplished  through  the  terms  in  the  primitive  equations 
which  have  been  presented  earlier. 


4.  Biological  Processes 

The  biological  model  is  a  system  of  N  coupled  partial  differential  equations  that  govern 
the  vertical  and  time  distribution  of  a  non-conservative  quantity,  and  is  of  the  form; 
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=  !Fb  +  S, 


(8) 


-—  +  {W+  Wbiology)  \^dz) 

where  5  is  a  non- conservative  quantity  (one  of  the  N  components  in  the  biological  model), 
w  is  the  vertical  advective  velocity,  wnology  is  the  vertical  sinking  rate  of  the  biological 
component,  Tb  is  the  cHmatological  forcing  term  and  S  is  the  source  or  sink  term  for  B. 

5.  Optical  Model  i  j  j  x- 

The  bio-optical  portion  of  the  model  simulates  the  wavelength-,  depth-  and  time- 

dependent  downwelhng  irradiance  field  between  300  and  4000  nm  with  a  spectral  resolu¬ 
tion  as  high  as  5  nm.  Both  the  direct  and  dilFuse  components  of  spectral  irradiance  at  the 
sea  surface  are  calculated  using  a  simple,  wavelength- dependent,  solar  irradiance  model  for 
direct-normal  and  diffuse  horizontal  irradiance  (Bird,  1984).  The  wavelength- dependent 
attenuation  of  the  subsurface  irradiance  field,  due  to  sea  water,  phytoplankton  and  dis¬ 
solved  organic  matter  is  simulated  using  a  modified  version  of  the  Sathyendranath  and 
Platt  (1986)  subsurface  solar  irradiance  model.  The  gradient  of  the  flux  of  the  subsur¬ 
face  solar  irradiance  is  incorporated  as  a  depth-dependent  energy  flux  which  balances  the 
phytoplankton  energy  uptake  and  the  thermal  energy  flux  (AT)  into  the  water  column. 

6.  Model  Implementation 

A  spectral  collocation,  method  was  used  to  solve  the  system  of  equations  that  define 
the  physical-bio-optical  model.  The  model  uses  a  leap-frog  timestepping  scheme  with 
occasional  leap-frog  trapazoidal  or  Adams-Bashforth  timestep  in  order  to  remove  the  effect 
of  timesphtting  from  the  solution.  This  provides  consistency  with  the  approach  used  to 
obtain  solutions  to  the  3-D  circulation  model  (Haidvogel  et  al.,  19916)  that  provided  input 
to  the  current  configuration  of  this  model  (Moisan,  1993).  A  more  complete  description  of 
this  spectral  technique  can  be  found  in  Haidvogel  et  al.  (1991a). 
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C.  Hood/Olson  PWP+NPZD  Coupled  Model 

Raleigh  Hood  and  Donald  Olson,  RSMAS,  University  of  Miami 

We  are  currently  working  with  a  coupled  PWP/NPZD  biophysical  model  (Price, 
Weller,  and  Pinkel  (1986)  mixed  layer  model  driving  a  nitrogen  (iV),  phytoplankton  (P), 
zooplankton  (Z)  and  detritus  {D)  ecosystem  model). 

The  four  compartment  biological  model  was  developed  by  Glen  Fherl,  Don  Olson  and 
Cabell  Davis  (unpubhshed).  The  basic  set  of  equations  are: 


dN 

dt 


((-^)gPZ  +  SdZ  +  eD  -  uIPN 


—  =  uIPN  -  gPZ  -  sP 
dt 


dt 


'ygPZ  -  dZ 


dt 


{l-C)gPZ  +  {l-6)dZ  +  sP  -  eD 


The  terms  are  defined  in  the  accompanying  table  along  with  a  schematic  box  diagram 
of  the  model  showing  the  pathways  and  the  directions  of  flows  of  nitrogen  between  the 

four  compartments. 

In  this  version  the  phytoplankton  growth  rate  is  Hnearly  dependent  upon  both  hght 
and  nutrients,  zooplankton  grazing  rate  is  hnearly  dependent  upon  phytoplankton  biomass, 
and  the  closure  is  the  simple  Z  form.  As  such,  equihbrium  solutions  are  easily  obtained. 
During  the  URIP  workshop  we  incorporated  Michaehs-Menten  saturation  functions  for 
production  versus  irradiance  and  production  versus  nitrogen  concentration,  and  a  modified 
Ivlev  zooplankton  grazing  function  (see  Physics  Working  Group  Report,  PWP+NPZD,  this 
volume).  We  also  plan  to  use  a  Z^  closure  in  future  simulations. 

The  essential  difference  between  NPZD  and  simpler  NPZ  models  (e.g.  Flierl  and 
Davis,  1993)  is,  obviously,  the  additional  D  compartment.  This  box  can  be  considered  a 
simple  “microbial  loop”,  where  the  detritus  reminerahzation  rate  (e)  is  the  rate  at  which 
bacteria  convert  detritus-nitrogen  to  dissolved-nitrogen.  Thus  in  the  NPZD  model  bac¬ 
terial  processes  are  exphcitly  represented.  Also,  the  D  compartment  allows  more  reahstic 

simulation  of  sinking. 

Our  NPZD  biological  model  runs  as  a  subroutine  in  a  modified  version  of  PWP. 
Nitrogen  concentrations  in  each  compartment  are  passively  mixed  (and  optionally  diffused) 
by  the  PWP  mixing  (and  diffusion)  subroutines.  As  with  temperature  and  salt,  the  mixing 
is  done  to  satisfy  the  stabihty  criterion  in  the  physical  model  (See  Price  et  al.,  1986).  Our 
version  of  PWP  has  been  modified  as  follows. 
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•  The  original  diffusion  subroutine  in  PWP  (forward  time  centered  space  or  FTCS)  has 
been  replaced  with  a  Crank-Nicholson  scheme  (Press  et  al.,  1986),  and  the  bottom 
boundary  condition  has  been  modified  to  allow  upward  diffusion  of  dissolved  nitrogen. 

•  An  additional  subroutine,  which  is  analogous  to  the  Crank-Nicholson  diffusion  sub¬ 
routine,  has  been  added  that  allows  detrital  sinking  and  loss  at  the  bottom  boundary. 

•  The  dual  wavelength  fight  attenuation/absorption  model  (Paulson  and  Simpson,  1977) 
has  been  modified  to  take  into  account  the  effects  of  phytoplankton  absorption.  The 
phytoplankton  absorb  only  blue  fight  according  to  a  constant  nitrogen-specific  diffuse 
attenuation  coefficient  (equivalent  to  a  Kc  —  0.014  w?  mgChla  ;  Kirk,  1983).  It 
is  assumed  that  only  10%  of  the  fight  absorbed  by  phytoplankton  is  converted  to 
photochemical  energy,  and  that  the  remaining  90%  ends  up  as  heat. 

•  We  have  incorporated  empirical  surface  forcing  formulae  from  Bleck  et  al.  (1989) 
that  provide  ideahzed  annual  cycles  in  insolation,  heat  loss  (latent,  sensible  and  back 
radiation)  and  wind  stress  for  the  North  Atlantic  ocean  as  functions  of  latitude. 

Coupled  Indian  Ocean  Model: 

In  collaboration  with  Jay  McCreary  at  NOVA  University  we  have  also  installed  &nNPZD 
ecosystem  model  in  a  2.5-layer  Indian  Ocean  circulation  model  (the  latter  is  described 
in  McCreary  et  al.,  1993).  A  surface  mixed  layer  embedded  in  the  upper  layer  of  the 
physical  model  determines  the  fight  (mixed  layer  depth)  and  nutrient  (entrainment  and 
detrainment)  environments  for  the  biological  model.  The  four  biological  compartments  are 
advected  horizontally  and  sinking  is  simulated  in  the  vertical. 
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N-P-Z-D  Model  Symbols  and  Units: 


Description 

Symbol 

Units 

Dissolved  nitrogen 

N 

Phytoplankton  nitrogen 

P 

fjM 

Zooplankton  nitrogen 

Z 

fjM 

Detritus  nitrogen 

D 

Light  intensity 

I 

dimensionless 

Phytoplankton  growth  rate 

u 

Zooplankton  Grazing  rate 

9 

Zooplankton  Death  rate 

d 

Detritus  remineralization  rate 

e 

Phytoplankton  scenscence  rate 

s 

Growth  efficiency  for  Z 

1 

dimensionless 

Assimililation  efficiency  for  Z 

c 

dimensionless 

Assimilation  efficiency  for  higher 
trophic  levelsf 

5 

dimensionless 

tgrowth  efficiency  for  higher  tronhic  levels  set  at  zero 
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